Enumeration of distinct mechanically stable disk packings in small systems 
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We create mechanically stable (MS) packings of bidisperse disks using an algorithm in which we 
successively grow or shrink soft repulsive disks followed by energy minimization until the overlaps 
are vanishingly small. We focus on small systems because this enables us to enumerate nearly all 
distinct MS packings. We measure the probability to obtain a MS packing at packing fraction 
4> and find several notable results. First, the probability is highly nonuniform. When averaged 
over narrow packing fraction intervals, the most probable MS packing occurs at the highest <j> and 
the probability decays exponentially with decreasing <j>. Even more striking, within each packing- 
fraction interval, the probability can vary by many orders of magnitude. By using two different 
packing-generation protocols, we show that these results are robust and the packing frequencies do 
not change qualitatively with different protocols. 

PACS numbers: 81.05.Rm, 82.70.-y, 83.80.Fg 



I. INTRODUCTION 



Inherent structures or potential energy minima are im- 
portant for determining the mechanical and dynamical 
properties of supercooled liquids and glasses p|. Re- 
cently, a host of computational studies have attempted to 
relate thermodynamic quantities in supercooled liquids 
to the number of inherent structures and vibrational mo- 
tions about them 0, 0, • It is often assumed in these 
calculations that all inherent structures at a given energy 
are equally probable. Similarly, statistical descriptions of 
granular media assume that all stable particle packings 
at a particular volume are equally likely 0, . Thus, it 
is important to examine under what conditions and to 
what extent that inherent structures in glassy systems 
and jammed packings of granular materials are equally 
probable. In this short article, we begin to address this 
question by enumerating nearly all mechanically stable 
(MS) packings in small 2d bidisperse systems. 

An innovative feature of this work and our other recent 
study is that we focus on systems containing small 
numbers of disks. We note that related studies of small 
hard disk systems have been carried out previously, but 
these have not investigated the MS packing probabilities 
0, 0. We confine our studies to small 2d systems for 
two key reasons. First, we are able to generate nearly 
all of the mechanically stable disk packings in these sys- 
tems. The number of MS disk packings grows exponen- 
tially with the number of particles N, but is finite for 
any finite N. In small systems, we are able to accurately 
determine the probability with which each MS packing 
occurs. In contrast, in large systems, only the most fre- 
quent MS packings are found. Second, we believe that 
understanding small jammed systems is crucial to devel- 
oping a theoretical explanation for slow stress and struc- 
tural relaxation in large glassy and amorphous systems 





FIG. 1: Snapshot of (a) a MS packing at <f> « 0.82 and (b) 
another MS packing (f> w 0.83 that is 10 6 times more frequent. 

|1C| . Our view of the relationship between jamming in 
small systems and glassy behavior in large systems will 
be more fully developed in a future publication [Tll | . 

We discuss several important results below. First, we 
find that the frequency with which mechanically stable 
packings occur is not uniform, instead it can vary by 
many orders of magnitude over the range of packing frac- 
tions where MS packings exist. This large variation in 
frequency occurs even though we do not target specific 
packings. Moreover, MS packings can have frequencies 
that differ by many orders of magnitude even over narrow 
packing-fraction intervals and there are no striking struc- 
tural differences between rare and frequent MS packings 
at similar (j) as shown in Fig. ^ Thus, it is not imme- 
diately obvious what variables set the frequency of MS 
packings. Also, we show below that even when we signif- 
icantly alter the protocol used to generate the MS pack- 
ings, the most frequent packings remain frequent and 
the rare ones remain rare. Thus, we suggest that al- 
though it is clear that the particular protocol chosen to 
generate the MS packings plays some role in determining 
the frequency distribution of MS packings |l2j. prominent 
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geometrical features of the configuration space can also 
strongly influence the frequency distribution. 

II. SIMULATION METHODS 
A. Generation of Mechanically Stable Packings 

We study mechanically stable packings of TV = 10 
disks that interact via the finite-range, pairwise additive, 
purely repulsive spring potential 



Vfa) = r(l - r«/ ff «) 2 e(<7«/r« - 1). (1) 



Here e is the characteristic energy scale, nj is the separa- 
tion between particles i and j, <Jij = (at + <Jj) /2 is their 
average diameter, and O(x) is the Heaviside step func- 
tion. We study 50-50 binary mixtures of particles with 
diameter ratio 1.4 to inhibit crystallization |7l ll5 . ll4|. ll5j. 
The particles have equal mass m and are confined to a 
square unit cell with periodic boundary conditions. We 
chose the small particle diameter a and a^/m/e as the 
reference length and time scales. 

Our focus is on configurations that are at a poten- 
tial energy minimum with infinitesimal overlaps. Since 
these configurations are in mechanical equilibrium and 
possess vanishingly small overlaps, we term them me- 
chanically stable (MS) packings. MS packings with no 
overlaps are equivalent to collectively jammed states 0] 
of hard disks. 

We employ a class of packing-generation protocols that 
involve successive compression or decompression steps 
followed by energy minimization 0, IT^ |. The process 
is initiated by choosing random initial positions for the 
particles at packing fraction <j> = 0.60, which is well be- 
low the minimum packing fraction at which MS pack- 
ings occur in 2d. The system is decompressed when the 
potential energy of the system at a local minimum is 
nonzero; otherwise, the system is compressed. The incre- 
ment by which the particle packing fraction is changed 
at each compression or decompression step is gradually 
decreased. After a sufficiently large number of steps a 
MS packing with infinitesimal overlaps is obtained. This 
process is performed for a large number of independent 
starting conditions to generate an ensemble of MS pack- 
ings. In this way, we can measure the probability to 
obtain a MS packing at a given <p. 

We employ two energy- minimization methods: (a) 
conjugate-gradient (CG) minimization algorithm or (b) 
molecular dynamics (MD) with dissipation proportional 
to local velocity differences. The conjugate-gradient 
method is a numerical scheme that begins at a given 
point in configuration space and moves the system to the 
nearest local potential e nerg y minimum without travers- 
ing any energy barriers |l8j . In contrast, molecular dy- 
namics with finite damping is not guaranteed to find the 
nearest local potential energy minimum since kinetic en- 
ergy is removed from the system at a finite rate. The 



system can thus surmount a sufficiently low energy bar- 
rier. In the molecular dynamics method, each particle i 
obeys Newton's equations of motion 
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where a, is the acceleration of particle i, 



(2) 

is the relative 



velocity of particles i and j, fy is the unit vector con- 
necting the centers of these particles, and b — 0.5 is the 
damping coefficient . Note that in our previous studies 
HQ HBI we used only the CG method. 



B. Classification of Mechanically Stable Packings 

In our numerical simulations we distinguish distinct 
mechanically stable disk packings by the lists of eigenval- 
ues of their dynamical matrices. For a pairwise additive, 
rotationally invariant potential Q the dynamical matrix 
is given by the expressions |l9j 
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and 



M; 



(4) 



where = dV/drij and Cij — d 2 V/dr?j. In the above 
relations the indices i and j refer to the particles, and 
a,/3 — x,y represent the Cartesian coordinates. For a 
system with Nf rattlers and N' — N — Nf particles form- 
ing a connected network the indices i and j range from 
1 to N', because the rattlers do not contribute to the 
potential energy. 

The dynamical matrix is symmetric and has dN' rows 
and columns, where d = 2 is the spatial dimension. 
Thus it has dN' real eigenvalues {rrii}, d of which are 
zero due to translational invariance of the system. In 
a MS disk packing, no set of particle displacements is 
possible without creating an overlapping configuration; 
therefore the dynamical matrix has exactly d(N' — 1) 
nonzero eigenvalues. In our simulations we use the crite- 
rion | mi | > m m i n = 10~ 6 for nonzero eigenvalues. 

We consider two MS packings to be the same if they 
have the same list of eigenvalues of the dynamical matrix. 
The eigenvalues are considered to be equal if they differ 
by less than the noise threshold m min for our calculations. 
Using the CG and MD methods, we have identified w 
1600 distinct MS packings for systems with 10 particles. 
(Packings with rattlers have been included in this count.) 

It is in general not true that each distinct MS pack- 
ing possesses a unique packing fraction <f>. However, we 
find that for these systems only at most a few percent 
of distinct MS packings share the same packing fraction. 
Thus, in the following we will associate a unique <fi with 
each MS packing to simplify the discussion. 
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FIG. 2: The probability to obtain a MS packing at packing 
fraction cj> for N = 10 particles using the MD method. The 
inset shows a magnified view of the probability on a logarith- 
mic scale over a narrow packing-fraction interval between 0.81 
and 0.82. 
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FIG. 3: The probabilities for the most frequent MS packings 
obtained from the CG packing-generation algorithm (filled 
squares) are compared to the probabilities for the same MS 
packings obtained from the MD method (filled circles). A 
similar comparison of less frequent MS packings obtained from 
the CG (open squares) and MD (open circles) methods are 
also shown. 



III. RESULTS 

The principal result of this work is that the frequency 
with which MS packings occur can vary dramatically — by 
many orders of magnitude from one packing to another. 
In Fig. |2 we show the probability P(<f>) — n((p)/N t to 
obtain a given MS packing at <f> for the MD energy min- 
imization method. n((p) is the number of MS packings 
obtained at cj> out of N t trials performed. For N = 10, we 
ran more than 30 x 10 6 trials for both energy minimiza- 
tion methods and obtained results that do not depend on 
N t . 

We notice several interesting features in Fig. [21 First, 
the probability distribution is not uniform. The most 
probable MS packings occur at large <j> « 0.83 and the 
least probable occur for tj> < 0.75. When averaged over 
narrow packing- fraction intervals, probability distribu- 
tion decays exponentially with decreasing <f> [jj. Even 
more striking, the probability is not monotonic in pack- 
ing fraction and is in fact noisy and difficult to predict. 
In the inset to Fig. |3 we find that even over a narrow 
range A(f>, the probability varies by more than five orders 
of magnitude, and this occurs over the entire range of </>. 

To understand the influence of the packing-generation 
protocol on our results, we have examined the proba- 
bilities for obtaining each MS packing using two energy 
minimization methods. We compared the 100 most fre- 
quent MS packings obtained from the CG method to the 
100 most frequent packings from the MD method. We 
find that ~ 80% of the MS packings were common to 
both sets; these are shown as filled symbols in Fig. [3] A 
similar comparison of probabilities for obtaining less fre- 
quent MS packings using the MD and CG methods is also 
displayed (using open symbols) in Fig. [3| We make sev- 



eral important observations. First, significant shuffling 
of frequent and rare packings does not occur — frequent 
packings remain frequent and rare ones remain rare. Sec- 
ond, within each probability grouping (i.e. filled versus 
open symbols) the more frequent MS packings in the CG 
method become more frequent in the MD method, while 
the less frequent packings from the CG method become 
more rare in the MD method. Although these results 
were obtained for small systems, we expect similar find- 
ings for large systems. 

The result that frequent packings become more fre- 
quent and rare ones become more rare when switching 
from CG to MD energy minimization can be explained 
in part by considering the rate at which energy is dis- 
sipated in the system. The MD energy minimization 
method utilizes a finite rate of energy dissipation. Thus, 
at various points during the packing generation process, 
the system can in principle possess enough kinetic energy 
to jump out of a shallow basin corresponding to a rare 
MS packing and into the basin of a more frequent pack- 
ing. This suggests that understanding the topography of 
configuration space surrounding MS packings is crucial 
to understanding the frequency of MS packings. 

Also, the result that frequent MS packings remain fre- 
quent and rare ones remain rare when switching from the 
MD to CG methods implies that geometrical features of 
configuration space strongly influence the frequency of 
MS packings. Correlations between the frequency with 
which MS packings occur and the shape and volume of 
basins near each MS packing is a topic of our current 
investigations [TT| . 
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IV. CONCLUSION 

We have enumerated nearly all of the distinct MS in 
small systems composed of 50-50 mixtures of bidisperse 
disks. This has allowed us to accurately measure the 
probability to obtain each MS packing. We find that the 
probability is not uniform — packings with large <fi are 
exponentially more likely than those at low <f>. More- 
over, even over narrow packing- fraction intervals, dis- 
tinct MS packings occur with frequencies that differ by 
many orders of magnitude. We have varied the packing- 
generation protocol and shown that these results do not 
change qualitatively, which suggests that geometrical fea- 
tures of configuration space strongly influence the pack- 
ing frequencies. These results suggest that further work 
should be carried out to explicitly test the Edwards' as- 



sumption that stable grain packings are equally probable 
and similar assumptions about inherent structures. An- 
other intriguing possibility is that only the most frequent 
MS packings are relevant for slow dynamics in glassy and 
jammed systems. However, even if this were true, we 
must first understand what determines MS packing fre- 
quencies in order to identify correctly the relevant set of 
frequent MS packings. 
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